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Abstract. We consider fractional order optimal control problems in which the 
dynamic control system involves integer and fractional order derivatives and 
the terminal time is free. Necessary conditions for a state/control/terminal- 
time triplet to be optimal are obtained. Situations with constraints present at 
the end time are also considered. Under appropriate assumptions, it is shown 
that the obtained necessary optimality conditions become sufficient. Numer- 
ical methods to solve the problems are presented, and some computational 
simulations are discussed in detail. 

1. Introduction. Fractional calculus generalizes the standard integral and differ- 
ential calculus to arbitrary order. If the order of the fractional derivative operator 
is an integer m, then we recover an m-fold integral when m is negative, and the 
classical derivative of order m when m is positive. Let x : [a,b] — > R be a func- 
tion, and let a be a real number representing the order of the integral/derivative 
of x. The most common fractional operators considered in the literature are the 
Riemann-Liouville fractional integrals (RLFI) and derivatives (RLFD), and the 
Caputo fractional derivatives (CFD): 

a I?x(t) = — - / (t - T) a - 1 x(r)dT, (left RLFI) 
"(a) Ja 

t I%x{t) = f (r - t) Q - 1 x(r)dr, (right RLFI) 
1 d n /"* 

a D?x(t) = y^^^ J a (* " r) n - a -^(r)dr, (left RLFD) 

tD?x{t) = ^-±- ) — J (r - t)»- a - l x{T)dT, (right RLFD) 



c a D?x(t) = / (t t)"-"-^ (r)dr, (left CFD) 

"(n - a) J a 

?D%x(t) = ] \ (r - tr^x^ (r)dr, (right CFD) 

T(n - a) J t 
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where n stands for [a\ +1 (n is the smallest integer larger than a), and T is the 
gamma function, that is, 



T(z) 



t z 



z > 0. 



If x and x^' , i = 1, . . . , n— 1, vanish at t = a, then a Dfx{t) — % Dfx{t), and if they 
vanish at t = b, then tD"x(t) — fD£x(t). The reader interested in the theory and 
applications of fractional calculus is referred to the books [11, 16]. For the state of 
the art on fractional order optimal control, see the recent book [15]. 

For numerical purposes, it is sometimes useful to approximate the fractional oper- 
ators as sums of integer order (standard) derivatives. The most common procedure 
is to replace the Riemann-Liouville fractional derivative by the series 



a D?x(t) = J2 



k=0 



a\ (t-a) k - a 
kjT(k + l-a) 



x (fe) (i), where 



(-l)*- 1 ar(fc-a) 
T(l - a)T(k + 1) 



(cf. [16]). This formula can be easily deduced by applying integration-by-parts to 
the integral 

1 d n rt 



x(t)cIt. 



r(n - a) dt n J a 

The main drawback is that to obtain good accuracy, we need higher-order deriva- 
tives of x, which restricts the class of admissible functions that one can consider. 
Recently, in [5], a new expansion formula was obtained, with the advantage that we 
only need the first derivative: 



a D?x(t) = A(a){t - a)' a x(t) + B(a)(t - a) 1 - a x{t) - C(a,p)(t - af-^Vpit), 

p=2 

where V p (t) is defined as the solution of the system 

iv p (t) = (i-p)(t- a y- 2 x(t), 

\V p (a) = 0, 

p = 2, 3, . . ., and the coefficients A, B and C are given by the formulas 



A(a) 



B(a) 



C(a,p) 



1 



r(l-a) 

1 

T(2 - a) 



1 + E 



p=2 



T(p - 1 + a) 
r(a)(p-l)! 



i 



r(p - 1 + a) 
^ r(a-l)p! _ 

r(p - 1 + a) 



r(2-a)r(a-l) 0-1)! 
For computational purposes, we truncate the sum and consider the finite expansion 

a D?x(t) ~ A(a, N){t - a)~ a x(t) + B(a, N){t - a) 1 - a x{t) 



Y,C{a,p)(t-af- 

p=2 



-p—a 



v P (t), (1) 
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where A(a, N) and B(a, N) are now denned by 

i N 

A(a, N) 



r(i 



1 



r(p - 1 + a) 



^r(a)(p-l)! 



B(a,N) 



1 



T(2 - a) 



E 



r(p-i 



=J r(a-l)p! 



We refer to [21] where expansion formulas with higher-order derivatives are ob- 
tained, and an error estimation is proven. In [20] analogous results are proven 
for the Riemann-Liouville fractional integral and in [19] for Hadamard fractional 
operators. For the right Riemann-Liouville fractional derivative, we have 

t D£x(t) ~ A(a, N)(b - t)- a x(t) - B{a, N)(b - tf^x^) 

N 

+ Y / C(a,p)(b-t) 1 -r~ a W p (t), (2) 



where W p is the solution of the differential equation 

(w p (t) = -(i-p)(b-ty- 2 x(t), 

\W p (b)=0. 

To approximate the Caputo fractional derivatives, we may use the formulas 



c a D?x{t) = a D?x(t) 



x^(a) 



and 



fc=0 



n-1 



r(fe-a + l) 



(t-a) 



k—a 



Cr>-> 



Dfx(t) = t D%x(t) - 



xW(b) 



(b - t) k - a 



fe=0 



r(fe-a+l) 



that relate Caputo and Riemann-Liouville fractional operators, and then use (1) 
and (2). In this article we are interested in investigating such ideas in the context 
of fractional order optimal control. An important tool is the integration by parts 
formula for Caputo fractional derivatives, which is stated in the following theorem. 

Theorem 1.1 (cf., e.g., [2]). Let a € (0,1), and x,y : [a, b] — > R be two functions 
of class C 1 . Then the following integration by parts formula holds: 



y{t) c a D?x{t)dt=[ t ll- a y(t)x{t)] b a + / x(t) t D?y(t)dt 



The text is organized as follows. In Section 2 we formulate the optimal con- 
trol problem under consideration and deduce necessary optimality conditions for 
it (Theorem 2.1). Another approach consists of using the approximation methods 
mentioned above, thereby converting the original problem into a classical optimal 
control problem that can be solved by standard computational techniques (Sec- 
tion 2.2). A generalization of the obtained results is considered in Section 3, where 
the lower bound of the cost integral is different from the lower bound of the frac- 
tional derivative (Theorem 3.1). Under some additional assumptions, the necessary 
optimality conditions are also sufficient (Theorem 4.1). We end with Section 5 
containing numerical computations. 
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2. Necessary optimality conditions. Let a € (0,1), a G K, L and / be two 

differentiable functions with domain [a, +oo) x R 2 , and <fi : [a,+oo)xR— s>Ra 
differentiable function. The fundamental problem is stated in the following way: 

minimize J(x,u,T) = / L(t,x(t),u(t))dt + <f>(T,x(T)) (3) 

J a 

subject to the control system 

Mx(t) + N c a D«x{t) = f(t, x(t), u(t)) (4) 
and the initial boundary condition 

x{a) = x a (5) 

with (M, N) ^ (0, 0) and x a a fixed real number. Our goal is to generalize previous 
works on fractional optimal control problems by considering the end time T free and 
the dynamic control system (4) involving integer and fractional order derivatives. 
For convenience, we consider the one-dimensional case. However, using similar 
techniques as the ones given here, the results can be easily extended to problems 
with multiple states and multiple controls. Later we consider the cases T and/or 
x(T) fixed. Here, T is a variable number with a < T < oo. Thus, we are interested 
not only on the optimal trajectory x and optimal control function u, but also on 
the corresponding time T for which the functional J attains its minimum value. We 
assume that the state variable x is differentiable and that the control u is piecewise 
continuous. When N — we obtain a classical optimal control problem; the case 
M = with fixed T has already been studied for different types of fractional order 
derivatives (see, e.g., [1, 3, 4, 8, 9, 22, 23]). In [10] a special type of the proposed 
problem is also studied for fixed T. 

Remark 1. In this paper the terminal time T is a free decision variable and, a 
priori, no constraints are imposed. For future research, one may wish to consider a 
class of fractional optimal control problems in which the terminal time is governed 
by a stopping condition. Such problems were recently investigated, within the 
classical (integer order) framework, in [13, 14]. 

2.1. Fractional necessary conditions. To deduce necessary optimality condi- 
tions that an optimal triplet {x, u, T) must satisfy, we use a Lagrange multiplier to 
adjoin the dynamic constraint (4) to the performance functional (3). To start, we 
define the Hamiltonian function H by 

H(t, x, u, A) = L(t, x, u) + \f(t, x, u), (6) 

where A is a Lagrange multiplier, so that we can rewrite the initial problem as 
minimizing 

J(x,u,T,X)= f [H(t,x,u,X)-X(t)[Mx(t)+N^Df L x(t)]] dt + 4>(T,x{T)), 

J a 

Next, we consider variations of the form 



x + <5.t, u + Su, T + ST, X + SX 
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with 8x(a) — by the imposed boundary condition (5). Using the well-known fact 
that the first variation of J must vanish when evaluated along a minimizer, we get 

f dM f)M f)M 

- X(t) (M5x(t) + N°D?Sx(t)j ^dt + ST [H(t, x, u, X) 

-X(t) (Mx{t) + N c a D?x(t))] t=T +^(T, x(T))6T+^(T, x(T)) (x(T)ST + Sx(T)) 

with the partial derivatives of H evaluated at (t, x(t), u(i), A(i)). Integration by 
parts gives the relations 



\{t)5x(t) dt = - Sx(t)X(t) dt + Sx(T)X(T) 



and 



X{t)^D?8x{t)dt = / 5x(t) t D%\{t)dt + 6x(T)[tI^ a \(t)}t=T. 



Thus, we deduce the following formula 
cT r 



5x 



f dH 
\~dx~ 



+ MX(t) - N t D^X{t) + 5 



dH 

du 



+5X 
8x(T) 



f dH 
[dX~ 



Mx(t)-N^ a D?x{t) 



dt 



MX(t)+N t I$T a X(t)-^(t,x(t)) 

d<t>, 



t=T 



+5T[ff(t,a:,u,A)-A(t)[Mx(t)+JV^i? t ol a:(t)]+^(t,a;(t))+^(t,a;(t))i(t)]t=r=0. 

Now, define the new variable 

Sx T = (x + 5x){T + 5T) - x{T). 

Because Sx(T) is arbitrary, in particular one can consider variation functions for 
which Sx(T) — 0. By Taylor's theorem, 

(x + 5x)(T + 8T) -(x + Sx)(T) = x(T)5T + 0(ST 2 ), 

where lim ^1 i s finite, and so Sx(T) — Sxt — x(T)ST + 0(ST 2 ). In conclusion, 
we arrive at the expression 



ST 



H(t,x,u,X) - NX{t)° D?x{t) + Nx{t) t I^T a X{t) + ^(t,x(t)) 



t=T 



dH 



. dH 
+Su— 
du 



Ox 
dt — 5xt 



6x (lE: + MA W - N tD%X(t) )+SX[^-~ Mx(t) - N°D?x(t) 



dH 



MX(t) + N t I±- a X(t) - ^(t,x(t)) 
dx 



+ 0{ST 2 ) = 0. 



t=T 



Since the variation functions were chosen arbitrarily, the following theorem is proven. 
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Theorem 2.1. If (x,u,T) is a minimizer of (3) under the dynamic constraint (4) 
and the boundary condition (5), then there exists a function X for which the triplet 
(x, u, A) satisfies: 

• the Hamiltonian system 



UlX(t) - N t D%X{t) = -™(t,x(t),u(t),\(t)) 
\Mx(t)+N c ;D?x(t) = %l(t,x(t),u(t),\(t)) 



(7) 



for all t € [a, T]; 
• the stationary condition 

dH 

du 



(t,x{t),u(t),X(t)) = 



(8) 



for all t G [a, T] ; 
• and the transversality conditions 

H(t,x(t),u(t), X(t)) ~N\{t) c a D?x{t) +Nx{t) t I 1 T - a \(t) + ^(t,x(t)) 



t=T 



MX(t) + N t I l T a X(t) - g 



where the Hamiltonian H is defined by (6). 



(t,x(t)) 



0; 



t=T 



(9) 



Remark 2. In standard optimal control, a free terminal time problem can be 
converted into a fixed final time problem by using the well-known transformation 
s = t/T (see Example 2). This transformation does not work in the fractional 
setting. Indeed, in standard optimal control, translating the problem from time t 
to a new time variable s is straightforward: the chain rule gives jjj — ^f^f- For 
Caputo or Riemann-Liouville fractional derivatives, the chain rule is not valid and 
such conversion is not possible. 

Some interesting special cases are obtained when restrictions are imposed on the 
end time T or on x(T). 

Corollary 1. Let (x,u) be a minimizer of (3) under the dynamic constraint (4) 
and the boundary condition (5). 

1. If T is fixed and x(T) is free, then Theorem 2.1 holds with the transversality 
conditions (9) replaced by 



MX{t) + N t I^r a X{t) - g^(t,x(t)) 



= 0. 



2. If x(T) is fixed and T is free, then Theorem 2.1 holds with the transversality 
conditions (9) replaced by 



H(t, x(t), u(t),X(t)) - NX(t)%D?x(t) + Nx(t)tI±~ a X(t) 



dcj) 
~dt 



(t,x(t)) 



0. 



t=T 



3. If T and x(T) are both fixed, then Theorem 2.1 holds with no transversality 
conditions. 

4- If the terminal point x(T) belongs to a fixed curve, i.e., x(T) — j(T) for 
some differentiable curve 7, then Theorem 2.1 holds with the transversality 
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conditions (9) replaced by 



H(t, x{t),u(t), A(i)) - N\(t) c a D?x(t) + Nx(t) t I^ a \(t) + -g{t, x{t)) 



j(t) (MX(t) + N t I*- a \(t) - ^(t, x(t)^j 



= 0. 



J t=T 



5. If T is fixed and x(T) > K for some fixed K g R, then Theorem 2.1 holds 
with the transversality conditions (9) replaced by 



MX{t) + N t I^r a X{t) - g^(t,x{t)) 



< o, 



(x(T) - K) 



M\(t) + N t I^ a X(t) - g(t, a;(t)) 



0. 



6. // x(T) is fixed and T < K for some fixed K g R, £/ien Theorem 2.1 holds 
with the transversality conditions (9) replaced by 



H(t, x(t), u(t),X(t)) - NX(t)°D?x(jt) + Nx(t) t I^- a \(t) + -£(t, x(t)) 



> 0, 



t=T 



H(t, x(t), u(t),X(t)) - N\{t) c a D«x{t) + Nx{t) t l]T a \{t) + -£(t, x(t)) 



dt 



x (T - K) = 0. 



Proof. The first three conditions are obvious. The fourth follows from 

Sx T = 7(T + ST) - 7 (T) = j(T)8T + 0{5T 2 ). 

To prove 5, observe that we have two possible cases. If x(T) > K, then Sxt may 
take negative and positive values, and so we get 



M\(t)+N t I±- a \(t) - g^{t,x(t)) 



On the other hand, if x(T) = K, then 8xt > and so 



M\(t) +N t I±- a \{t) - g^(t,x(t)) 



0. 



< 0. 



t=T 



The proof of the last condition is similar. 



□ 



Case 1 of Corollary 1 was proven in [8] for (M, N) = (0, 1) and = 0. Moreover, if 
a = 1, then we obtain the classical necessary optimality conditions for the standard 
optimal control problem (see, e.g., [7]): the Hamiltonian system 

(±(t) = 2g(t,x(t),u(t),\(t)), 
\\(t) = -™(t 1 x(t),u(t),\(t)); 

OH 

the stationary condition — — (i, x(t), u(t), X(t)) = 0; and the transversality condition 
A(T) =0. 
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2.2. Approximated integer order necessary optimality conditions. Using 
approximation (1) up to order K, we can transform the original problem (3)-(5) 
into the following classical problem: 

i-T 



minimize J(x,u,T) 



L(t, x(t), u(t)) dt + <j)(T, x(T)) 



subject to 



x(t) 



f(t, x(t), u(t)) - NA(t - a)- a x(t) + Ef =2 NC p (t - a)*-P-°%(t) 
M + NB(t — a) 1_Q 



and 



V p (t) = (1 - p)(t - af- 2 x{t), p = 2,...,K, 
\x(a) = x a , 

\V p (a)=0, p = 2,...,K, 



(10) 



where A = A(a,K), B = B(a,K) and C v — C(a,p) are the coefficients in the 
approximation (f ). Now that we are dealing with an integer order problem, we can 
follow a classical procedure (see, e.g., [12]), by defining the Hamiltonian H by 



H = L(t,x,u) + 



Ai (/(t, x, u) - NA(t - a)- a x + £f =2 NC p (t - af-^V^ 



M + NB(t - a) 



l-a 



K 



+ ^A p (l- J3 )(t-a)P- 2 x. 

Let A = (Ai, A 2 , . . . , Xk) and x = (x, V2, ■ ■ ■ , Vk)- The necessary optimality condi- 
tions 

dH 



dH 



du 



= 0, 



- dH 

result in a two point boundary value problem. Assume that (T* ,x* ,u*) is the 
optimal triplet. In addition to the boundary conditions (10), the transversality 
conditions imply 



d: 



-(T*,x*(T)) 



<5x<r 



H(T*,x*(T),u*(T),\*(T)) + ^(T*,x*(T)) 



6T = 0, 



where tr denotes the transpose. Because V p , p = 2, . . . ,K, are auxiliary variables 
whose values V P (T), at the final time T, are free, we have 



Ap(T) 



t=T 



0, p = 2,. 



The value of Ai(T) is determined from the value of x(T). If x(T) is free, then 

dx 



Ai(T) = ^-\t=T- Whenever the final time is free, a transversality condition of the 



fori 







completes the required set of boundary conditions. 
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3. A generalization. The aim is now to consider a generalization of the optimal 
control problem (3)-(5) studied in Section 2. Observe that the initial point t = a 
is in fact the initial point for two different operators: for the integral in (3) and, 
secondly, for the left Caputo fractional derivative given by the dynamic constraint 
(4). We now consider the case where the lower bound of the integral of J is greater 
than the lower bound of the fractional derivative. The problem is stated as follows: 

minimize J(x,u,T) = [ L(t, x(t),u(t)) dt + <j>(T 7 x(T)) (11) 



under the constraints 

Mx{t) + N°D*x{t) = f(t, x(t), u(t)) and x{A) = x A , (12) 
where (M, N) ^ (0, 0), xa is a fixed real, and a < A. 

Remark 3. We have chosen to consider the initial condition on the initial time A 
of the cost integral, but the case of initial condition x(a) instead of x(A) can be 
studied using similar arguments. Our choice seems the most natural: the interval 
of interest is [A, T] but the fractional derivative is a nonlocal operator and has 
"memory" that goes to the past of the interval [A, T] under consideration. 

Remark 4. In the theory of fractional differential equations, the initial condition is 
given at t = a. To the best of our knowledge there is no general theory about unicity 
of solutions for problems like (12), where the fractional derivative involves x(t) for 
a < t < A and the initial condition is given at t = A. Unicity of solution is, however, 
possible. Consider, for example, QD"x(t) = t 2 . Applying the fractional integral 
to both sides of equality we get x(t) = x(0) + 2t 2+a /T(3 + a) so, knowing a value 
for x(t), not necessarily at t — 0, one can determine x(0) and by doing so x(t). 
A different approach than the one considered here is to provide an initialization 
function for t G [a, A]. This initial memory approach was studied for fractional 
continuous-time linear control systems in [17] and [18], respectively for Caputo and 
Riemann-Liouville derivatives. 

The method to obtain the required necessary optimality conditions follows the 
same procedure as the one discussed before. The first variation gives 



J A 





I A 



^—5x + ^—Su + ^-SX - SX (Mx{t) + N c a D?x(t)) 
ox ou oX 



X(t) (M8x(t) + N^D^Sx(t) 



dt + ^(T, x(Tj) {x{T)5T + Sx(Tj) 
ox 



+ -£(T, x{T))8T + ST [H(t, x, u, X) - X(t) (Mx(t) + N c a D?x{t))] t=T . 
where the Hamiltonian H is as in (6). Now, if we integrate by parts, we get 

T pT 

X(t)Sx(t) dt = - 5x(t)X(t) dt + 5x(T)X(T) 
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and 



\(t%D?8x(t)dt= \(t%D?5x{t)dt- X{t) c a D?5x{t)dt 



5x(t) t D%\(t) dt + [Sx(t) t I^r a X(t)]iZa 



Sx(t) t D%X(t) dt 



[6x(t) t I^X(t)Y=t 

A 

8x{t)\ t D«X{t) - t D a A X(t)]dt 



Sx{t) t D%X(t) dt 

+ 5x(T)[ t I^r a X(t)] t=T - 8x(a)[ a I^ a X(a) - J^ a X{a)]. 
Substituting these relations into the hrst variation of J, we conclude that 



dH 



dx 



dH 



— — I- MA - N t D^X Sx + —Su + — - Mx - N°D?x 5X 



du 



dH 



dX 



dt 



N 



-(t,x)] 



t=T 



f A 

\ Sx[ t D%X- t D1X}dt- Sx[MX + N t I^r a X- , 

J a ® X 

+ ST[H{t, x, u, A) - A [Mi + N°D?x] + ^-{t, x) + ^-(t, x)x] t=T 

dt ox 

+ NSx(a)[ a I^ a X(a) - a ljf a A(a)] = 0. 

Repeating the calculations as before, we prove the following optimality conditions. 

Theorem 3.1. If the triplet (x,u,T) is an optimal solution to problem (11)— (12), 
then there exists a function X for which the following conditions hold: 

• the Hamiltonian system 

d H 

MX(t) - N t D%X(t) = —g^(t, *(*). u{t),X{t)) 
Mx{t)+N c a D?x(t) = —(t,x(t),u(t),X(t)) 

for all t G [A, T], and t D%X(t) - t D%X(t) = for all t £ [a, A]; 

• the stationary condition 

OH 

— (t,x(t),u{t),X{t))=0 

for all t e [A, T]; 

• the transvcrsality conditions 



H(t, x(t), u(t),X(t)) - NX(t)^Dfx(t) + Nx(t) t I^r a X(t) + -£(t, x(t)) 



MX(t)+N t I^r a X(t)-^(t,x(t)) 
ox 



dt 

= 0, 



= 0, 



t=T 



t=T 



[ f 4- Q A(t)- t /i-«A(t)] t=a = 0; 
with the Hamiltonian H given by (6). 

Remark 5. If the admissible functions take fixed values at both t = a and t = A, 
then we only obtain the two transversality conditions evaluated at t = T. 
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4. Sufficient optimality conditions. In this section we show that, under some 
extra hypotheses, the obtained necessary optimality conditions are also sufficient. 

Theorem 4.1. Let (x,u,X) be a triplet satisfying conditions (7)-(9) of Theo- 
rem 2.1. Moreover, assume that 

1. L and f are convex on x and u, and <f> is convex in x; 

2. T is fixed; 

3. X(t) > for all t £ [a,T] or f is linear in x and u. 
Then (x,u) is an optimal solution to problem (3) -(5). 

Proof. From (7) we deduce that 



dL — — — df 

— (t,x(t),u{t)) = -MX(t) + N t D$X(t) - X(t)-^-(t,x(t)Mt))- 



Using (8), 



^(t,x(t),u(t)) = -A(t)|£(t,2(t),5(t)) 



and (9) gives [MX(t) + N t L? a X(t) - §f (*, z(i))]t=T = 0. Let (x,u) be admissible, 
i.e., let (4) and(5) be satisfied for (x, u). In this case, 

J(x, u) — J(x, u) 



[L(t, x(t), u(t)) - L(t, x(t),u(t))} dt + 0(T, x(T)) - <j>(T, x(T)) 



> 



^(t,x(t)Mt))(x(t) -x(t)) + ^(t,x(t),u(t))(u(t) -u(t)) 



dt 



+ -^(T,x(T))(x(T)-x(T)) 



-MX(t)(x(t) - x(t)) + N t D%X(t)(x(t) - x{t)) 



X(t)-^-(t,x(t),u(t))(x(t) -x{t)) - X{t)-±(t,x{t),u{t)){u{t) ~u{t)) 



dx 



th 



dt 



+ -£(T,x(T))(x(T)-x(T)). 



Integrating by parts, and noting that x(a) = x(a), we obtain 



J(x, u) — J(x, u) 



X(t) 



M (x(t) - 3?(t)) + N { c a D«x{t) - G a D?x{t)) 
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|£ (t,x(t),u(t)) (x(t)-x(t)) - |£ (t,x(t),u(t)) (u(t)-u(t)) 



dt 



■£(t,x(t)) - MX(t) - N t I^ a \(t) 



(x(T)-x(T)) 



t=T 



A(i) [f(t, x(t), u(t)) - / (t,x(t),u(t))} - X(t)^(t,x(t),u(t)) (x(t) - x(t)) 



■\(t)^(t,x(t),u(t))Ht)-u(t)) 



dt 



rp 

> I A(t) 



^ (t,5(t),5(t)) (as(t) - z(t)) + ^- (f,S(t),H(*)) (u(t) - S(t)) 



df df 

-± {t,x(t),u(t)) (x(t)-x(t)) - ^ (t, *(*),«(*)) («(*) -u(0) 



dt 



0. 



□ 



Remark 6. If the functions in Theorem 4.1 are strictly convex instead of convex, 
then the minimizer is unique. 

5. Numerical treatment and examples. Here we apply the necessary condi- 
tions of Section 2 to solve some test problems. Solving an optimal control problem, 
analytically, is an optimistic goal and is impossible except for simple cases. There- 
fore, we apply numerical and computational methods to solve our problems. In 
each case we try to solve the problem either by applying fractional necessary con- 
ditions or by approximating the problem by a classical one and then solving the 
approximate problem. 



5.1. Fixed final time. We first solve a simple problem with fixed final time. In 
this case the exact solution, i.e., the optimal control and the corresponding optimal 
trajectory, is known, and hence we can compare it with the approximations obtained 
by our numerical method. 

Example 1. Consider the following optimal control problem: 

J(x, u) = [ (tu(t) - (a + 1)x(t)f dt — ► min 
Jo 

subject to the control system 

x(t) + £ D?x{t) =u(t)+t 2 
and the boundary conditions 

2 



x(0) = 0, x(l) = 



r(3 + a)' 



The solution is given by 
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because J(x,u) > for all pairs (x,u) and x(0) = 0, x(l) — fpq^g ; x(t) = u(t) 

and ^Dfx(t) = t 2 with J(x,u) = 0. It is trivial to check that (x, u) satisfies the 
fractional necessary optimality conditions given by Theorem 2.1 /Corollary 1. 

Let us apply the fractional necessary conditions to the above problem. The 
Hamiltonian is H = (tu — (a + 2)x) 2 + \u + At 2 . The stationary condition (8) 
implies that for t ^ 



and hence 



Finally, (7) gives 



u[t) = — - 



A 2 . > -r - , _> , 



H = -— + — — arA + rA, t + 0. (13) 



| x(t) + $D?x{t) = + *&x(t) + t 2 , \ x(0) = 0, 
\-A(i)+ t D?A(t) = s±2A(t), \ x (l) = r^j. 

At this point, we encounter a fractional boundary value problem that needs to be 
solved in order to reach the optimal solution. A handful of methods can be found 
in the literature to solve this problem. Nevertheless, we use approximations (1) and 
(2), up to order N, that have been introduced in [5] and used in [10, 21]. With 
our choice of approximation, the fractional problem is transformed into a classical 
(integer order) boundary value problem: 

'*(*) = [(H 1 - At~<*) x(t) + e^ =2 c p t^-v p (t) - ^ + 1 2 ] TT5 ^ 

V p (t) = (l-p)tP- 2 x(t), p=2,...,N 

a(*) = [(a(i ty<* sf) X (t) - Ef= 2 c p (i - ty-p-w p (t)] 1+B[1 \ t) ^ 

[W p (t) = -(1 - p)(l - t)f- 2 A(i), p — 2, . . . , N 
subject to the boundary conditions 

'as(0)=0, a:(l) = r^5J' 
V p (0) = 0, p = 2,...,N, 
W P (1) = 0, p = 2,...,N. 

The solutions are depicted in Figure 1 for N — 2, N = 3 and a = 1/2. Since the 
exact solution for this problem is known, for each N we compute the approximation 
error by using the maximum norm. Assume that x(U) are the approximated values 
on the discrete time horizon a = to,ti, . . . , t n . Then the error is given by 

E = max(|x(i;) -x(U)\). 

i 

Another approach is to approximate the original problem by using (1) for the 
fractional derivative. Following the procedure discussed in Section 2, the problem 
of Example 1 is approximated by 

J(x,u) = / (tu — (a + 2)x) 2 dt — > min 
Jo 

subject to the control system 

\x(t)[l + B(a, N)^-*] + A(a, N)t- a x(t) - Y^Li C(a,p)t 1 ~P- a V p (t) = u(t) + t 2 
\V p (t) = (l-p)tP- 2 x(t) 
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(c) x(t),N = 3 (d) u(t),N = 3 

Figure 1. Exact solution (solid lines) for the problem in Exam- 
ple 1 with a = 1/2 versus numerical solutions (dashed lines) ob- 
tained using approximations (1) and (2) up to order N in the frac- 
tional necessary optimality conditions. 



and boundary conditions 

as(0)=0, x(l)= 2 , V p (0)=0, p = 2,3,...,N. 
1 [6 + a) 

The Hamiltonian system for this classical optimal control problem is 

Xij-Aja, N)t- a x + Y,p=2 Cj^p^-P^Vp + u + t 2 ) 



H={tu-{a + 2)xf 



1 + B(a,N)t 1 - a 

N 

^(l-p)^- 2 A p x. 



Using the stationary condition ^ = 0, we have 



, , a + 2 , , Ai (t) 

u(t) = ——x(t) ^ , r-i for t + 0. 

w t w 2t 2 (l + B(a,N)t 1 ~ a ) r 
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Finally, the Hamiltonian becomes 

N N 

H = o \ 2 1 + cl> 1 x\ 1 +J2 ( f > pVpX 1 +<f> N+1 \ 1 +J2( 1 -p) tP ~ 2x \' MO, (14) 



p=2 



p=2 



where 



and 



4&(l + B{a,N)t 1 -°'y 



a + 2 - A(a, N)^-" 
t(l + B(a,N)t 1 - a ) ' 



l + S(a,A^)^ 



(15) 



(16) 



The Hamiltonian system x 



dH \ 



ff , gives 



i(t) = 20 o (i)Ai(i) + 0i(i)«(t) + Ep=2^(W(<) + 
y p = (1 - P )tP- 2 x(t), p = 2,...,N 

\\ = (*) + E^= 2 (P - l)^ 2 A P 

.A p = -^(t)A 1 (t), p = 2,...,iV 



subject to the boundary conditions 



x(0) = 

V p (0) = 0, p = 2,. 



. iV 



- r(3+ Q ) 
A P (1) = 0, p = 2,...,iV. 



This two-point boundary value problem was solved using Matlab's bvp4c built-in 
function for N = 2 and N = 3. The results are depicted in Figure 2. 

5.2. Free final time. The two numerical methods discussed in Section 5.1 are now 
employed to solve a fractional order optimal control problem with free final time T. 

Example 2. Find an optimal triplet (x(-),u(-),T) that minimizes 

rT 



J(x,u)= {tu - (a + 2)x) 2 dt 
Jo 

subject to the control system 

x(t) + oD?x{t) =u(t)+t 2 

and boundary conditions 

x(0)=0, x(T) = l. 

An exact solution to this problem is not known and we apply the two numerical 
procedures already used with respect to the fixed final time problem in Example 1. 

We begin by using the fractional necessary optimality conditions that, after ap- 
proximating the fractional terms, result in 

f x(t) = [(s±2 _ At-") x(t) + Ef =2 C p ^-"V p (t) - # + t 2 
V p (t) = (l-p)tP- 2 x(t), p=2,...,N 

X(t) = [(A(l - - s±2) X(t) - Ef =2 C p (l - ty-r>-*W p {t) 
[W p (t) = -(1 - p)(l - t)P- 2 X(t), p = 2, . . . , N 



l + Bt 1 - 



i+B(i~ty- 



1(5 
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Figure 2. Exact solution (solid lines) for the problem in Exam- 
ple 1 with a = 1/2 versus numerical solutions (dashed lines) ob- 
tained by approximating the fractional order optimal control prob- 
lem using (1) up to order N and then solving the classical necessary 
optimality conditions with Matlab's bvp4c built-in function. 



subject to the boundary conditions 

x(0)=0, x(T) = l, 
F p (0)=0, p = 2,...,N, 
W p (T) = 0, p = 2,...,N. 

The only difference here with respect to Example 1 is that there is an extra unknown, 
the terminal time T. The boundary condition for this new unknown is chosen 
appropriately from the transversality conditions discussed in Corollary 1, i.e., 

[H(t, x, u, A) - \{t)CD?x{t) + x{t) t I^ a \{t)] t=T = 0, 

where H is given as in (13). Since we require A to be continuous, t/^ _Q A(i)|t = T = 
(cf. [16, pag. 46]) and so A(T) = 0. One possible way to proceed consists in 
translating the problem into the interval [0, 1] by the change of variable t = Ts 
[6]. In this setting, either we add T to the problem as a new state variable with 
dynamics T(s) — 0, or we treat it as a parameter. We use the latter, to get the 



FRACTIONAL ORDER OPTIMAL CONTROL PROBLEMS 



17 




(a) x(t),N = 2 



(b) u(t),N ■■ 



Figure 3. Numerical solutions to the free final time problem of 
Example 2 with a = 1/2, using fractional necessary optimality 
conditions (dashed lines) and approximation of the problem to an 
integer order optimal control problem (dash-dotted lines). 



following parametric boundary value problem: 



x{s) 



[(^-A(Tsr a )x(s)+T," =2 CATs) 1 - p - a VAs)-^$p + (Ts) 2 ]T 



l + B(Ts) 1 - c 

V p (s) = T(l - p)(Ts)P' 2 x(s), p = 2,...,N, 

w x _ [(A(l-T s )"°-^)A( s )-Ef =2 C P (l-T s ) 1 -"- c 



W p (s)]T 



l+B^-Ts) 1 - 

[W p (s) = -T(l-p)(l-T s )P- 2 X( s ), 



subject to the boundary conditions 



P = 2,...,N, 




'ac(O) = 0, 

V p (0) = 0, p = 2,...,N, 
.W P (1) = 0, p = 2,...,N, 

This parametric boundary value problem is solved for N = 2 and a = 0.5 with 
Matlab's bvp4c function. The result is shown in Figure 3 (dashed lines). 

We also solve Example 2 with a = 1/2 by directly transforming it into an integer 
order optimal control problem with free final time. As is well known in the classical 
theory of optimal control, the Hamiltonian must vanish at the terminal point when 
the final time is free, i.e., one has H\ t= T = with H given by (14) [12]. For N = 2, 
the necessary optimality conditions give the following two point boundary value 
problem: 

'x(t) = 20 o (*)Ai(i) + 0i(i)a(t) + <h(t)V*(t) + fa(t) 
V 2 = -x(t) 

X 1 = - ( j) l (t)X 1 {t) + x(t) 
.A 2 = -&(f)Ai(t), 

where 4>o(t) and <fii(t) are given by (15) and and <j>$(f) by (16) with p = N = 2. 
The trajectory x and corresponding u are shown in Figure 3 (dash-dotted lines). 
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